Plane geometry and convexity 
of polynomial stability regions-^ 



Didier Henrion^'^, Michael Sebek^ 
February 2, 2008 



Abstract 

The set of controllers stabilizing a linear system is generally non-convex in the 
parameter space. In the case of two-parameter controller design (e.g. PI control or 
static output feedback with one input and two outputs), we observe however that 
quite often for benchmark problem instances, the set of stabilizing controllers seems 
to be convex. In this note we use elementary techniques from real algebraic geometry 
(resultants and Bezoutian matrices) to explain this phenomenon. As a byproduct, 
we derive a convex linear matrix inequality (LMI) formulation of two-parameter 
fixed-order controller design problem, when possible. 
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1 Introduction 

Despite its elementary formulation, the problem of fixed-order controller design for linear 
time-invariant systems remains mostly open. Especially scarce are numerically efficient 
computer-aided control system design algorithms in the fixed-order case, sharply con- 
trasting with the large number of tools available to solve static state feedback design or 
dynamical output feedback design with controllers of the same order as the plant. Math- 
ematically, fixed-order controller design can be formulated as a non-convex non-smooth 
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optimization problem in the parameter space. To the best of our knowledge, randomized 
algorithms are amongst the most efficient numerical methods to cope with this class of 
difficult problems. See [H [13] for computer experiments supporting this claim, thanks to 
public-domain Matlab packages (HIFOO and the Randomized Control System Toolbox). 

This note was motivated by the observation, made by the first author during a workshop 
at AIM in August 2005 [8j, that 6 out of the 7 two-dimensional instances of static output 
feedback (SOF) design problems found in the database COMPleib seem to be convex. 
Further motivation was provided by the excellent historical survey [7] on D-decomposition 
techniques, previously studied in deep detail in [T2] and [1]. In [7] the authors describe 
the intricate geometry of two-dimensional stability regions with the help of illustrative 
examples. Quite often, the stability regions represented in these references seem to be 
convex. 

In this note, we use basic results from real algebraic geometry to detect convexity of the 
stability region in the two-parameter case (including PI controllers, PID with constant 
gain, SOF design with one input two outputs or two inputs one output). We also derive, 
when possible, a linear matrix inequality (LMI) formulation of the stability region. 



2 Problem statement 



We consider a parametrized polynomial 

p{s,k) = po{s) + kipi{s) + k2P2{s) (1) 

where the Pi{s) G M[s] are given polynomials of s G C and the /cj G M are parameters. We 
assume, without loss of generality, that the ratio Pi{s) / p2{s) is not a constant. 

Define the stability region 

S = {keR^ : p{s, k) stable} 

where stabihty is meant in the continuous-time sense, i.e. all the roots of p(s, k) must lie 
in the open left half-plan^. 

We are interested in the following problems: 

• Is stability region S convex ? 

• If it is convex, give an LMI representation 

S = {keR^ : Fo + Fiki + ^2^2 ^ 0} 

when possible, where the Fi are real symmetric matrices to be found, and >- means 
positive definite. 



^Similar results can be derived for discrete-time unit disk stability or other semialgebraic stability 
domains of the complex plane, but this is not covered here. 
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2.1 Example: PI controller design 



Let b{s)/a{s) with a{s),b{s) G R[s] denote the transfer function of an open loop plant, 
and consider a proportional integral (PI) controller ki/s + k2 in a standard negative 
feedback configuration. The closed-loop characteristic polynomial ([1]) is then p{s, k) = 
sa{s) + {ki + k2s)h{s) hence po{s) = sa{s), pi{s) = b{s) and P2{s) = sb{s). 



2.2 Example: static output feedback 



Given matrices A E M"^", B e M"""", C G M^''", we want to find a matrix K G W'^^'p 
such that the closed-loop matrix A + BKC is stable. When mp = 2, the characteristic 
polynomial ([T]) writes p{s, k) = det(s/„ — A — BKC) hence po{s) = det(sJ„ — A), pi{s) = 
det(sJ„ - B[ki, 0]C) and p2(s) = det(s/„ - B[0, k2]C). 



3 Hermite matrix 



The Routh-Hurwitz criterion for stability of polynomials has a symmetric version called 
the Hermite criterion. A polynomial is stable if and only if its Hermite matrix, quadratic 
in the polynomial coefficients, is positive definite. In control systems terminology, the 
Hermite matrix is a particular choice of a Lyapunov matrix certifying stability [TT]. Al- 
gebraically, the Hermite matrix can be defined via the Bezoutian, a symmetric form of 
the resultant [HI Section 5.1.2]. 

Let a{u), b{u) be two polynomials of degree n of the indeterminate u. Define the Bezoutian 
matrix Bu{a, b) as the symmetric matrix of size n with entries bij satisfying the linear 
equations 

a{u)b{v) — a{v)b{u) 



EE" 



u 

i=i j=i 



The polynomial r„(a, b) = det Bu{a, b) is the resultant of a{u) and b{u) with respect to u. 
It is obtained by eliminating u from the system of equations a{u) = b{u) = 0. 



The Hermite matrix of p{s) is defined as the Bezoutian matrix of the real part and the 
imaginary part of p{juj): 

Pr{uj'^) = Rep{juj) 
Lvpiiuj"^) = lmp{juj) 

that is, H{p) = B^{j)ji{(jj'^),upj{(jj'^)). Let us assume that p{s) is monic, with unit leading 
coefficient. The Hermite stability criterion can be formulated as follows. 



Lemma 1 Polynomial p{s) is stable if and only if H{p) >- 0. 



Proof: The proof of this result can be found in for example. It can also be proved 
via Cauchy indices and Hermite quadratic forms for counting real roots of polynomials, 
see [21 Section 9.3] □. 

By construction, the Hermite matrix of parametrized polynomial ([1]) 
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H{p{s,k)) = H{k)= J2 H.n^,k\%^ 

jl,«2=0 

is quadratic in k. Therefore, the Hermite criterion yields a quadratic matrix inequality 
formulation of the stability region: 

S = {keR^ : H{k) y 0}. 

Quadratic matrix inequalities, a generalization of bilinear matrix inequalities, typically 
generate non-convex regions. For example, the scalar quadratic inequality kf — 1 > 
models a disconnected, hence non-convex set. 

Surprisingly, it turns out that S, even though modeled by a quadratic matrix inequahty, 
is often a convex set for practical problem instances. Here are some examples. 



3.1 Examples: static output feedback 

Consider the 7 two-parameter SOF problems found in the database COMPleib [LQl, la- 
belled AC4, ACT, AC17, NNl, NN5, NN17 and HEl. Stability regions are represented as shaded 
gray areas on Figures [T] to [3 Visual inspection reveals that 6 out of 7 stability regions 
seem to be convex. The only apparently nonconvex example is HEl. 




Figure 1: AC4. Figure 2: AC7. Figure 3: AC 17. 

In the remainder of the paper we will explain why such planar stability regions are likely 
to be convex, and how we can constructively derive their LMI formulations when possible. 



A 



Figure 4: NNl. 



Figure 5: NN5. 



Figure 6: NN17. 




4 Rational boundary of the stability region 

Define the curve 

C = {A; G : p{juj, k) = 0,w eR} 

which is the set of parameters k for which polynomial k) has a root along the boundary 
of the stability region, namely the imaginary axis. Studying this curve is the key idea 
behind the D-decomposition approach [7]. The curve partitions the plane (^1,^2) into 
regions in which the number of stable roots of p{s, k) remains constant. The union of 
regions for which this number is equal to the degree of p{s, k) is the stability region S. 
Hence the boundary of S is included in curve C. 

Note that p{juj, k) = for some w G M if and only if 

Pniuj^, k) = poRiuj^) + kipmiuj^) + A;2P2i?,(^^) = 
ujpiiuj'^^k) = upQi{uj) + kiujpuiuj) + k2Ujp2i{uj) = 



Recall that we denote by r^{qi,q2) the resultant of polynomials qi{uj) ,q2{uj) obtained by 
eliminating the scalar indeterminate u. From the definition of the Hermite matrix, it 
holds 

Hk) = r^(pi^(cu^ k),ujpi{iu^ k)) = det H{k) (2) 



from which the imphcit algebraic description 

C = {k: h{k) = 0} 

follows. 

Lemma 2 The determinant of the Hermite matrix can he factored as 

h{k) = mg{kf 

where l{k) is affine, and g{k) is a generically irreducible polynomial. 

Proof: The result follows from basic properties of resultants: h{k) = r^{pj^{uj'^,k),ujpi{uj'^,k)) = 
ruj{,PR{^^^,k),uj'^)r^{pR{uj'^,k),pi{uj'^,k)) = r^{pR{u,k),u)r^{pR{LJ,k),pi{uj,k)Y. Take g{k) = 
ru{pR{uj, k),pi{u, k)). Since Pr{uj, k) is affine in k it follows that l{k) = r^^pR^u, k),Lj) is 
affine in /c.D 

The curve can therefore be decomposed as the union of a line and a simpler algebraic 
curve 

C = Cug = {k: l{k) = 0}U{k: g{k) = 0}. 
The equation of line C was already given in the proof of Lemma [2], namely 

Kk) = r^{pR{uj),uj) = pr{0, k) = Por{0) + kipiR^O) + A;2P2fi(0). 

The defining polynomial of the other curve component Q can be obtained via the formula 

gik) = r^{pR{uj, k),pi{uj, k)). 

From the relations 

_ u:pii{u'^) ujp2i{uj'^) 
we derive a rational parametrization of Q: 

' k,{uj') 
_ k^iuj') 

which is well-defined since by assumption Pi{s)/p2{s) is not a constant. From this 
parametrization we can derive a symmetric linear determinantal form of the implicit 
equation of this curve. 
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Lemma 3 The symmetric affine pencil 

G{k) = B^{qi, ga) + kiB^{q2, Qo) + k2B^{qi, Qq). 
is such that Q = {k : detG'(A;) = 0}. 



Proof: Rewrite the system of equations (I3l) as 

a{uj'^,k) = gi(u;^) — fcigo(^^^) = 

and use the Bezoutian resultant to ehminate indeterminate u and obtain conditions for a 
point (fci, ^2) to belong to the curve. The Bezoutian matrix is B^{a, h) = B^{qi — kiqQ, — 
k2qo) = B^{qi,q2) + kiB^{q2,qo) + ^2-8(91, go)- Linearity in k follows from bilinearity of 
the Bezoutian and the common factor go-D 

Finally, let C{k) = diag {l{k), G{k)} so that curve C can be described as a determinantal 
locus 

C = {k : det Cik) = 0}. 



5 LMI formulation 

Curve C partitions the plane into several connected components, that we denote by Si for 
t = l,...,N. 

Lemma 4 If C{k) >- for some point k in the interior of Si for some i then Si = {k : 
C{k) ^ 0} is a convex LMI region. 

Proof: Follows readily from the afiine dependence of C{k) on k and from the fact that 
the boundary of Si is included in C. □ 

Convex sets which admit an LMI representation are called rigidly convex in [9]. Rigid 
convexity is stronger than convexity. It may happen that Si is convex for some z, yet 
C{k) is not positive definite for points k within Si. 

Lemma 5 Stability region S is the union of sets Si containing points k such that H{k) >- 
0. 

Proof: Follows readily from Lemma [Tjn 

Note that it may happen that Si is convex LMI for some i, yet H{k) is not positive definite 
for points k within Si. 

Corollary 1 If H{k) >- and C{k) )~ Q for some point k in the interior of Si, then Si is 
an LMI region included in the stability region. 
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Proof: Combine Lemmas H] and [5l □ 

Quite often, on practical instances, we observe that S = Si for some i is a convex LMI 
region. 

Practically speaking, once curve C is expressed as a determinantal locus, the search of 
points k such that C{k) y can be formulated as an eigenvalue problem, but this is out 
of the scope of this paper. 



6 Examples 
6.1 Example 1 



-8 -6 -4 -2 2 4 6 8 

Figure 8: Vishnegradsky's degree 3 polynomial. Rational curve C with convex LMI sta- 
bility region S (in gray) 



4 



2 



-2 



-4 



As mentioned in [7j, Vishnegradsky in 1876 considered the polynomial p{s,k) = + 
kis'^ + k2S + 1 and concluded that its stability region S = {k : ki > 0, kik2 > 1} is convex 
hyperbolic. 



The Hermite matrix of p{s, k) is given by 

H{k)- 

and hence after a row and column permutation the quadratic matrix inequality formula- 



's 1 

kik2 - 1 

1 ki 
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tion 



S = {k: 



k2 1 
1 ki 
kik2 - 



explains why the region is convex. Indeed, the determinant of the 2-by-2 upper matrix, 
affine in k, is equal to the remaining diagonal entry, which is here redundant. The stability 
region can therefore be modeled as the LMI 



S = {k: 



ki 1 
1 k2 



see Figure [HI 



6.2 Example 2 

100 




Figure 9: NNl SOF example. Rational curve C with convex LMI stability region S (in 
gray) 



Consider problem NNl from fiO\. for which pn(s) = — 13), = s{s — 5), p2{s) = s + 1 



in ([ID. We have PRiuj"^) = — fcio;^ + k2 and P/(u;^) 

k2(-13-5ki + k2) 



13 — 5ki + k2 and 



H{k) 









-k2 



ki{-13 - 5ki + k2) 




-k2 


ki 



Hence h{k) = det H{k) = k2{-l3ki - /c2 - 5kf + fcifc2)^ and then l{k) = /c2, g{k) = 
— 13/ci — ^2 — 5kl + kik2. A rational parametrization of the curve Q = {k : g{k) = 0} is 
given by 

k2{uj^) = u\uj^ + 13)/{uj^-5) 
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from which we derive the symmetric affine determinantal representation Q = {k : det G{k) 
0} with 

' 169 + 65ki - 18k2 13 + 5ki " 
13 + 5ki 1 - ki 



Gik) 



The pencil representing C is therefore 

C{k) = 



k2 

169 + 65fci - I8k2 13 + 5ki 
13 + 5^1 1 - ki 



We can check that S = {k : C{k) 0} is a convex LMI formulation of the stability 
region represented on Figure HI Compare with Figure [9] where we represent also the curve 
C = {k: det C{k) =0}. 



6.3 Example 3 




-2' ' ' ' — ' ' ' ' ' ' ' ' 1 

-4.5 -4 -3.5 -3 -2.5 -2 -1.5 -1 -0.5 0.5 1 

k, 



Figure 10: Francis' example. Rational curve C with convex LMI stability region S (in 
gray) 



This example, originally from Francis (1987), is also described in [7]. A SISO plant 
(s — l)(s — 2)/(s + l)(s^ + s + 1) must be stabilized with a PI controller ki + k2/s. 
Equivalently, po{s) = s{s + l)(s^ + s + 1), pi{s) = s{s — l)(s — 2) and P2{s) = (s — l)(s — 2) 
in (P. 



For these values we obtain 

2k2 





Cik) 







14 + 28A;i - 5Ak2 
-20 - 40A;i + 18A;2 
2 + Aki 





-20 - 40A;i + 18A;2 
77 - 53ki + 36^2 
-11 + 5A;i 





2 + 4A;i 
-ll + 5fci 
5 + ki 
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and the LMI stability region S = {k : C{k) >~ 0} represented on Figure fTOl together with 
the quartic curve C = {k : det C{k) = 0}. 



6.4 Example 4 

Consider P Example 14.4] for which po{s) = + 2s^ + lOs^ + 10s + 14 + 2a, pi{s) = 
2s^ + 2s — 3/10, P2{s) = 2s + 1, with a G M a parameter. 




Figure 11: Ackermann's example with a = 1. Rational curve C and stability region S 
consisting of two disconnected regions (in gray), the one including the origin being LMI. 



We obtain C{k) = diag {l{k), G{k)} with 
l{k) = 140 + 20a - 3ki + 10^2 



G{k) 



7920 + 4860a + 400a2 
+ (-1609 - 60a)A;i + (-270 + 200a)A;2 
-8350 - 2000a + 1430fci + 130^2 
900 + 200a - 130A;i 



8370 - 1230A;i - 100A;2 * 
-900 + lOOfci 100 



where symmetric entries are denoted by stars. 



When a = 1, the stability region consists of two disconnected components. The one 
containing the origin /ci = /c2 = is the LMI region {k : C{k) >- 0}, see Figure [TTl 

When a = 0, the stability region S is the non-convex region represented on Figure \T2[ 
The LMI region {k : C{k) >- 0} is not included in S in this case. 



1 1 



-10 -5 5 10 15 20 

Figure 12: Ackermann's example with a = 0. Rational curve C and non-convex stability 
region S (in gray). 

7 Conclusion 

In this paper we have explained why the planar stability region of a polynomial may be 
convex with an explicit LMI representation. This is an instance of hidden convexity of a set 
which is otherwise described by intersecting (generally non-convex) Routh-Hurwitz minors 
sublevel sets or by enforcing positive definiteness of a (generally non-convex) quadratic 
Hermite matrix. 

Practically speaking, optimizing a closed-loop performance criterion over an LMI formula- 
tion of the stability region is much simpler than optimizing over the non-linear formulation 
stemming from the Routh-Hurwitz minors or the Hermite quadratic matrix inequality. 

Convexity in the parameter space was already exploited in [HI [2] in the context of PID 
controller design. It was shown that when the proportional gain is fixed, the set of integral 
and derivative gains is a union of a finite number of polytopes. 

Extension of these ideas to the case of more than 2 parameters seems to be difficult. 
The problem of finding a symmetric affine determinantal representation of rationally 
parametrized surfaces or hypersurfaces is not yet well understood, to the best of our 
knowledge. For example, in the simplest third degree case p{s, k) = + kis^ + + fcs, 
how could we find four symmetric real matrices Aq, Ai, A2, satisfying det(y4o + Ai/ci + 
742^2 + ^43/^3) = kik^ - ^2 ? 
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